Set STAN Environment

install.packages("cmdstanr", repos = c("https://mc-stan.org/r-packages/", getOption("repos")))

## The two are used for STAN result analysis.
install.packages("bayesplot")
install.packages("posterior")
library(cmdstanr)
library(bayesplot)
library(posterior)

## for bayesplot plotting
color_scheme_set("brewer-Spectral")

## If you have cmdstan, let cmdstanr know it.
set_cmdstan_path(path = paste(Sys.getenv("HOME"),
  "softwares",
  "cmdstan-2.23.0", sep = "/"))

## set the parameter for STAN sampling with multi-core support.
options(mc.cores = 3)

Example 1: Multivariate normal with missing data

Data

Our observations are from \(\mathcal{N}(0, \Sigma)\), where \(\Sigma\) is unknown.

\[\begin{pmatrix} 1 & -1 & -1 & 1 & 2 & 2 & -2 & -2 & * & * & * & *\\ 1 & -1 & 1 & -1 & * & * & * & * & 2 & 2 & -2 & -2 \end{pmatrix}\]

Each column corresponds to one sample, and \(*\) are missing data.

How we infer \(\Sigma\) using STAN

  • load the model
camel_mdl <- cmdstan_model("camel/camel2.stan", compile = TRUE)
## Model executable is up to date!
camel_mdl$print()
## /*
##  * Camel: Multivariate normal with structured missing data
##  * http://www.openbugs.net/Examples/Camel.html
##  *
##  */
## transformed data {
##   vector[2] Y[4];
##   real Y1[4];  // missing y2
##   real Y2[4];  // msising y1
## 
##   matrix[2, 2] S;
## 
##   S[1, 1] = 1;
##   S[1, 2] = 0;
##   S[2, 1] = 0;
##   S[2, 2] = 1;
## 
##   Y[1, 1] = 1.;
##   Y[1, 2] = 1.;
##   Y[2, 1] = 1.;
##   Y[2, 2] = -1.;
##   Y[3, 1] = -1.;
##   Y[3, 2] = 1.;
##   Y[4, 1] = -1.;
##   Y[4, 2] = -1.;
## 
##   Y1[1] = 2.;
##   Y1[2] = 2.;
##   Y1[3] = -2.;
##   Y1[4] = -2.;
## 
##   Y2[1] = 2.;
##   Y2[2] = 2.;
##   Y2[3] = -2.;
##   Y2[4] = -2.;
##  vector[2] mu = rep_vector(0.0, 2);
## }
## 
## 
## parameters {
##  // vector[2] mu;
##   cov_matrix[2] Sigma;
## }
## 
## model {
##  // mu ~ normal(0, 10);
##  // Sigma ~ inv_wishart(3, S);
##   for (n in 1:4) {
##      Y[n] ~ multi_normal(mu, Sigma);
##  }
##   Y1 ~ normal(mu[1], sqrt(Sigma[1, 1]));
##   Y2 ~ normal(mu[2], sqrt(Sigma[2, 2]));
##  // Use Jeffery's prior: p(Sigma) proportiol to |Sigma|^{-3/2}
##  target += -1.5 * log(determinant(Sigma));
## }
## 
## generated quantities {
##  real<lower=-1,upper= 1> rho;
##   rho = Sigma[1, 2] / sqrt(Sigma[1, 1] * Sigma[2, 2]);
## 
##  real y2condy1[4];
##  real y1condy2[4];
##  real muy2condy1[4];
##  real muy1condy2[4];
## 
##  for (i in 1:4) {
##      muy2condy1[i] = mu[2] + sqrt(Sigma[2,2]) * rho * (Y1[i] - mu[1])/sqrt(Sigma[1,1]);
##      y2condy1[i] = normal_rng(muy2condy1[i],
##                                               sqrt((1-rho*rho) * Sigma[2,2]));
##  }
## 
##  for (i in 1:4) {
##      muy1condy2[i] = mu[1] + sqrt(Sigma[1, 1]) * rho * (Y2[i]- mu[2]) / sqrt(Sigma[2,2]);
##      y1condy2[i] = normal_rng(muy1condy2[i],
##                                                       sqrt((1-rho*rho) * Sigma[1,1]));
##  }
## }
  • sampling the parameter
camel_sampler <- camel_mdl$sample(data = NULL,
  seed = 355113,
  parallel_chains = 4,
  show_message = TRUE,
  ## increase adapt delta, default is 0.80
  adapt_delta = 0.85,
  iter_sampling = 1000
  ## output_dir = "./"
)
## Running MCMC with 4 parallel chains...
## 
## Chain 1 Iteration:    1 / 2000 [  0%]  (Warmup) 
## Chain 1 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 1 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 1 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 1 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 1 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 1 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 2 Iteration:    1 / 2000 [  0%]  (Warmup) 
## Chain 2 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 2 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 2 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 2 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 2 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 3 Iteration:    1 / 2000 [  0%]  (Warmup) 
## Chain 3 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 3 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 3 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 3 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 3 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 3 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 3 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 3 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 3 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 3 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 3 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 4 Iteration:    1 / 2000 [  0%]  (Warmup) 
## Chain 4 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 4 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 4 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 4 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 4 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 4 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 4 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 4 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 4 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 4 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 4 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 4 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 1 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 1 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 1 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 1 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 1 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 1 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 1 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 1 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 1 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 1 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 1 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 1 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 2 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 2 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 2 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 2 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 2 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 2 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 2 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 2 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 2 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 2 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 2 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 2 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 2 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 3 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 3 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 3 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 3 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 3 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 3 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 4 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 4 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 4 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 4 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 4 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 1 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 1 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 1 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 2 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 2 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 2 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 3 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 3 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 3 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 3 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 4 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 4 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 4 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 4 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 1 finished in 0.5 seconds.
## Chain 2 finished in 0.5 seconds.
## Chain 3 finished in 0.5 seconds.
## Chain 4 finished in 0.5 seconds.
## 
## All 4 chains finished successfully.
## Mean chain execution time: 0.5 seconds.
## Total execution time: 0.6 seconds.
camel_draws <- camel_sampler$draws()
  • show the posterior results
bayesplot::mcmc_combo(camel_draws,
  pars = c("rho"), widths = c(1, 3),
  combo = c("dens_overlay", "trace")
 )

Example 2: Finite normal mixture model

load STAN script for both simulation and inference

mixture_model <- cmdstan_model("mixture/normal_mixture_k_order.stan",
  compile = TRUE
)
## Model executable is up to date!

Generate the simulation data

N <- 1200
N1 <- 400
N2 <- 400
musig <- 10.0
alpha <- c(1.0, 1.0, 1.0)
sigma <- c(1.0, 1.0, 1.0)

sim_sampler <- mixture_model$sample(
  data = list(
    K = 3,
    N = N,
    N1 = N1,
    N2 = N2,
    y = rep(1, N),
    musig = musig,
    alpha = alpha,
    sigma = sigma
  ),
  seed = 355113,
  chains = 1,
  iter_sampling = 1,
  parallel_chains = getOption("mc.cores", 1),
  show_message = TRUE,
  init = list(
    list(
      theta = c(N1 / N, N2 / N, 1 - N1 / N - N2 / N),
      mu = c(-5.0, 0.0, 5)
      ## sigma = c(0.5, 1, 1)
    )
  ),
  fixed_param = TRUE
)
## Running MCMC with 1 chain...
## 
## Chain 1 Iteration: 1 / 1 [100%]  (Sampling) 
## Chain 1 finished in 0.1 seconds.
sim_draws <- sim_sampler$draws()
## dimnames(sim_draws)$variable[1:10]

gy <- c(sim_draws[1, 1, grep("gy", dimnames(sim_draws)$variable)])
hist(gy, breaks = 50)

Sampling with STAN

sim_sampler <- mixture_model$sample(
  data = list(
    K = 3,
    N = N,
    N1 = N1,
    N2 = N2,
    y = gy,
    musig = musig,
    alpha = alpha,
    sigma = sigma
  ),
  seed = 355113,
  chains = 4,
  iter_sampling = 1000,
  fixed_param = FALSE,
  ## init = list(
  ##   list(
  ##     theta = c(0.4, 0.3, 0.3),
  ##     mu = c(-4.0, 0.0, 4),
  ##     sigma = c(10, 10, 10)
  ##   ),
  ##   list(
  ##     theta = c(0.2, 0.3, 0.5),
  ##     mu = c(-10.0, 0.0, 10),
  ##     sigma = c(10.0, 10, 10.0)
  ##   )
  ## ),
  parallel_chains = getOption("mc.cores", 1),
  ## output_dir = "./",
  show_message = TRUE
)
## Running MCMC with 4 chains, at most 3 in parallel...
## 
## Chain 1 Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 1 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 1
## Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 1 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 1
## Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 1 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 1
## Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 1 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 1
## Chain 2 Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 2 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 2 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 2 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 2 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 2
## Chain 3 Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 3 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 3 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 3 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 3 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 3
## Chain 3 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 3 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 3 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 3 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 3
## Chain 2 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 1 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 3 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 2 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 3 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 1 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 2 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 3 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 1 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 2 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 3 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 1 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 2 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 3 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 1 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 2 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 3 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 1 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 2 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 3 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 1 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 2 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 3 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 1 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 2 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 3 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 1 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 2 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 2 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 3 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 3 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 1 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 1 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 2 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 1 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 3 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 2 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 1 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 3 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 2 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 1 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 2 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 3 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 1 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 2 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 3 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 1 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 2 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 3 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 1 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 2 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 3 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 1 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 2 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 3 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 1 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 2 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 3 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 1 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 2 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 2 finished in 9.7 seconds.
## Chain 4 Iteration:    1 / 2000 [  0%]  (Warmup)
## Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 4 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 4
## Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 4 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 4
## Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 4 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 4
## Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
## Chain 4 Exception: normal_lpdf: Location parameter is inf, but must be finite! (in '/var/folders/jv/thfh42cj1zv2k4csp7sk7qg40000gn/T/Rtmp20RcUG/model-175fb62371ccd.stan', line 28, column 3 to column 49)
## Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
## Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
## Chain 4
## Chain 3 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 1 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 1 finished in 10.3 seconds.
## Chain 3 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 3 finished in 10.6 seconds.
## Chain 4 Iteration:  100 / 2000 [  5%]  (Warmup) 
## Chain 4 Iteration:  200 / 2000 [ 10%]  (Warmup) 
## Chain 4 Iteration:  300 / 2000 [ 15%]  (Warmup) 
## Chain 4 Iteration:  400 / 2000 [ 20%]  (Warmup) 
## Chain 4 Iteration:  500 / 2000 [ 25%]  (Warmup) 
## Chain 4 Iteration:  600 / 2000 [ 30%]  (Warmup) 
## Chain 4 Iteration:  700 / 2000 [ 35%]  (Warmup) 
## Chain 4 Iteration:  800 / 2000 [ 40%]  (Warmup) 
## Chain 4 Iteration:  900 / 2000 [ 45%]  (Warmup) 
## Chain 4 Iteration: 1000 / 2000 [ 50%]  (Warmup) 
## Chain 4 Iteration: 1001 / 2000 [ 50%]  (Sampling) 
## Chain 4 Iteration: 1100 / 2000 [ 55%]  (Sampling) 
## Chain 4 Iteration: 1200 / 2000 [ 60%]  (Sampling) 
## Chain 4 Iteration: 1300 / 2000 [ 65%]  (Sampling) 
## Chain 4 Iteration: 1400 / 2000 [ 70%]  (Sampling) 
## Chain 4 Iteration: 1500 / 2000 [ 75%]  (Sampling) 
## Chain 4 Iteration: 1600 / 2000 [ 80%]  (Sampling) 
## Chain 4 Iteration: 1700 / 2000 [ 85%]  (Sampling) 
## Chain 4 Iteration: 1800 / 2000 [ 90%]  (Sampling) 
## Chain 4 Iteration: 1900 / 2000 [ 95%]  (Sampling) 
## Chain 4 Iteration: 2000 / 2000 [100%]  (Sampling) 
## Chain 4 finished in 10.2 seconds.
## 
## All 4 chains finished successfully.
## Mean chain execution time: 10.2 seconds.
## Total execution time: 20.0 seconds.

View the results.

ey <- sim_sampler$draws()
bayesplot::mcmc_combo(ey[, , grep("mu", dimnames(ey)$variable)],
  widths = c(1, 3),
  combo = c("dens_overlay", "trace"))

If we set unordered \(mu\) in STAN

uy <- umixture_sampler$draws()
bayesplot::mcmc_combo(uy[, , grep("mu", dimnames(uy)$variable)],
                      widths = c(1, 3),
                      combo = c("dens_overlay", "trace"))